/*
 * DevCalc.cpp
 *
 *  Created on: 14.11.2012
 *      Author: Frantisek Schneider
 */

#include "GenCalc.h"
#include "CsvLoader.h"
#include <boost/random.hpp>
#include "tbb/parallel_for.h"
#include <stdlib.h>
#include "tbb/task_scheduler_init.h"
#include "ParallelWorker.h"

using namespace std;

GenCalc::GenCalc(data* data, Configuration* conf)
{
	this->datas = data;
	this->bestFittnes = WORSE_FITTNES;
	this->sumSelection = 0;
	this->mutationRate = MUTATION_RATE;
	this->threadCount = conf->threadCount;
	this->bestMember = NULL;
	this->conf = conf;
}

GenCalc::~GenCalc()
{
	for(unsigned int i=0; i<population.size(); i++)
	{
		delete[] population.at(i).gens;
	}

	delete[] bestMember;
	population.clear();
}

void GenCalc::rankPopulation()
{
	sumSelection = 0;
	double maxSelection = 100;				// the best one on the top gets 100% chance to be selected

	for (unsigned int i = 0; i < population.size(); i++)
	{
		population.at(i).selectionFitness = maxSelection;
		sumSelection += population.at(i).selectionFitness;
		maxSelection /= DECREASE_RATE;
	}
}

void GenCalc::applyConfigurationState(double* &gens)
{

	if(conf->p != -1)
	{
		gens[0] = conf->p;
	}

	if (conf->c != -1)
	{
		gens[2] = conf->c;
	}
}

double* GenCalc::getBestMember()
{
	return this->bestMember;
}

void GenCalc::initPopulation()
{
	boost::mt19937 rng;

	boost::uniform_real<float> pRange(0, 1);
	boost::uniform_real<float> cgRange(-1, 0);
	boost::uniform_real<float> cRange(0, 1);
	boost::uniform_real<float> dtRange(conf->dtmin, conf->dtmax);
	boost::uniform_real<float> kRange(conf->kmin, conf->kmax);

	boost::variate_generator<boost::mt19937&, boost::uniform_real<float> > gen1(rng, pRange);
	boost::variate_generator<boost::mt19937&, boost::uniform_real<float> > gen2(rng, cgRange);
	boost::variate_generator<boost::mt19937&, boost::uniform_real<float> > gen3(rng, cRange);
	boost::variate_generator<boost::mt19937&, boost::uniform_real<float> > gen4(rng, dtRange);

	boost::variate_generator<boost::mt19937&, boost::uniform_real<float> > gen5(rng, kRange);

	member newMember;

	for (int i = 0; i < conf->populationCount; i++)
	{
		double* gens = new double[5];

		gens[0] = gen1();		// p
		gens[1] = gen2();		// cg
		gens[2] = gen3();		// c
		gens[3] = gen4();		// dt
		gens[4] = gen5();		// k

		applyConfigurationState(gens);

		newMember.gens = gens;
		population.push_back(newMember);
	}
}


struct FitnessComparator
{
    inline bool operator() (const member& struct1, const member& struct2)
    {
        return (struct1.fitness < struct2.fitness);
    }
};

int GenCalc::pickMember()
{
	double value = sumSelection * (double) rand() / RAND_MAX;
	double sum = 0;


	for (unsigned int i = 0; i < population.size(); i++)
	{
		sum += population.at(i).selectionFitness;
		if (sum >= value)
		{
			return i;
		}
	}
	return rand() % population.size();
}

void GenCalc::crossover(double* result, double* first, double* second)
{
	for (int i = 0; i < 5; i++)
	{
		double p = (double) rand() / RAND_MAX;
			if (p < CROSS_RATE)
			{
				result[i] = (first[i] + second[i]) / 2;
			}
			else
			{
				result[i] = first[i];
			}
	}
}

void GenCalc::mutate(double* gens)
{
	for (int i = 0; i < 5; i++)
	{
		double ppb = (double) rand() / RAND_MAX;
		if (ppb < MUTATION_RATE)
		{
			gens[i] = gens[i] + ppb;
		}
	}
}
double GenCalc::getBestFitness()
{
	return this->bestFittnes;
}


void GenCalc::compute()
{
	for (int i = 1; i < conf->generationCount + 1; i++)
	{
		nextGeneration();
		nextEvolution(i);
	}
}

void GenCalc::nextGeneration()
{
		member* newPopulation = new member[population.size()];

		int size = population.size() /2 ;

		for (int i = 0; i < size; i++)
		{
			int member1 = pickMember();
			int member2;
			while (member1 == (member2 = pickMember()));

			newPopulation[i*2].gens = new double[5];
			newPopulation[i*2 + 1].gens = new double[5];

			crossover(newPopulation[i*2].gens, population.at(member1).gens, population.at(member2).gens);
			crossover(newPopulation[i*2 +1].gens, population.at(member2).gens, population.at(member1).gens);

			mutate(newPopulation[i*2].gens);
			mutate(newPopulation[i*2 + 1].gens);

			applyConfigurationState(newPopulation[i*2].gens);
			applyConfigurationState(newPopulation[i*2+1].gens);
		}


	for(unsigned int i=0; i<population.size(); i++)
	{
		delete[] population.at(i).gens;
	}

	population.clear();

	for(int i=0; i<size; i++)
	{
		population.push_back(newPopulation[i*2]);
		population.push_back(newPopulation[i*2+1]);
	}

	delete[] newPopulation;
}

void GenCalc::writeFittnesAtPosition(double fitness, int index)
{
	population.at(index).fitness = fitness;
}

void GenCalc::nextEvolution(int generation)
{

//	tbb::parallel_for(tbb::blocked_range<int>(0, population.size()),			// g++ na hydre neumi lambdu :(
//			[=] (const tbb::blocked_range<int>& r)
//			{
//				for (int i = r.begin(); i < r.end(); i++)
//				{
//					fitness(population.at(i));
//				}
//			});
//	tbb::task_scheduler_init terminate();


//	for (unsigned int i = 0; i < population.size(); i++)
//	{
//		population.at(i).fitness = fitness(population.at(i));
//	}

	tbb::task_scheduler_init init(threadCount);
	ParallelWorker pw(population, this, datas);
	tbb::parallel_for(tbb::blocked_range<int>(0, population.size()), pw);
	tbb::task_scheduler_init terminate();


	sort(population.begin(), population.end(), FitnessComparator());

	rankPopulation();

	if (population.at(0).fitness < this->bestFittnes)
	{
		if(this->bestMember != NULL)
		{
			delete[] this->bestMember;
		}

		this->bestMember = new double[5];
		this->bestMember[0] = population.at(0).gens[0];
		this->bestMember[1] = population.at(0).gens[1];
		this->bestMember[2] = population.at(0).gens[2];
		this->bestMember[3] = population.at(0).gens[3];
		this->bestMember[4] = population.at(0).gens[4];			//memcpy nefacha kdyz je to moc rychly :/

		this->bestFittnes = population.at(0).fitness;

		printf("[%d] : Best values (fitt=%f):  p=%f, cg=%f, c=%f , dt=%f, k=%f \n", generation,
				this->bestFittnes,
				bestMember[0],
				bestMember[1],
				bestMember[2],
				bestMember[3],
				bestMember[4]);
	}
}


bool GenCalc::checkTimeCondition(double k, double dt, double ecg)
{
	double minutes = k * dt * ecg * MINUTE;
	if (-5 < minutes || minutes > 5)
	{
		return false;
	}
	return true;
}

double GenCalc::fitness(member inputMember)
{
	double res = 0;
	double* gens = inputMember.gens;

	double p = gens[0];
	double cg =gens[1];
	double c = gens[2];
	double dt =gens[3];
	double k = gens[4];

	for (int i = 0; i < datas->length; i++)
	{
		double ecg = datas->items[i].ecg;
		if (checkTimeCondition(k, dt, ecg))
		{
			res = WORSE_FITTNES;
			return res;
		}

		double left_side = p * datas->items[i].bg+ cg * datas->items[i].bg * ecg + c;
		double index_time = datas->items[i].time + dt + k * dt * ecg;
		double ig;

		if (this->getIg(index_time * (24 * 60 * 60), &ig) == 1)
		{
			double right_side = ig;
			res += fabs(right_side - left_side);
		}
	}

	return res;
}

int GenCalc::getIg(double time, double *ig)
{
	*ig = 0;

	int t = (int) time;
	int timeStart = this->find_time(t);

	if (timeStart == -1 || timeStart >= datas->length - 1)
	{
		return 0;	 //FALSE
	}

	*ig = datas->items[timeStart].ig + (datas->items[timeStart + 1].ig -  datas->items[timeStart].ig) * (time - t);

	return 1;
}

int GenCalc::find_time(unsigned int time)
{

	int start, end;
	int midIndex = 0;			// binary search
	start = 0;
	end = datas->length;

	while (start <= end)
	{
		midIndex = (start + end) / 2;
		if (datas->items[midIndex].t == time)
		{
			return midIndex;
		}
		else if (datas->items[midIndex].t < time)
		{
			start = midIndex + 1;
		}
		else
		{
			end = midIndex - 1;
		}
	}
	return -1;
}
